Snotel
AnnaClaire Marley
6/18/2020
This RMarkdown conducts initial time series analyses on two SNOTEL sites in the Chuska Mountains located in the Navajo Nation. It then compares SNOTEL data to remote sensing data such as SNODAS and PRISM to understand which products best account for snow.
Read in packages
First explore NAs
## # A tibble: 7 x 3
## variable n_miss pct_miss
## <chr> <int> <dbl>
## 1 snow_depth_mm 44 1.09
## 2 air_temp_max_C 29 0.718
## 3 air_temp_av_C 25 0.619
## 4 air_temp_min_C 25 0.619
## 5 swe_mm 15 0.371
## 6 precip_mm 1 0.0248
## 7 Date 0 0
The horizontal black bars in the lower left indicate the number of NAs total for each variable, with the variable name shown to the right of the corresponding horizontal bar
The vertical black lines with black dots indicate the variables between which the frequency of co-occurring NAs is indicated by the vertical black bars
Daily basic trends
## Warning: Removed 1 rows containing missing values (geom_path).## Warning: Removed 7 rows containing missing values (geom_path).
ENSO
## Warning: Removed 1 rows containing missing values (geom_path).
blue = la nina years
red = el nino years
## Picking joint bandwidth of 25.1
Monthly basic trends
Water year trends
Compare to SNODAS data (extracted for the same lat/long as snotel site)
Daily data
## Joining, by = "Date"## Picking joint bandwidth of 19.6
## Picking joint bandwidth of 69.3## Warning: Removed 5 rows containing non-finite values (stat_density_ridges).
## Warning: Removed 2049 rows containing missing values (geom_path).
## Adding missing grouping variables: `month`
## Warning: Removed 2055 rows containing missing values (geom_path).
## Warning: Removed 5 rows containing non-finite values (stat_smooth).## Warning: Removed 5 rows containing missing values (geom_point).
## Warning: Removed 5 rows containing non-finite values (stat_smooth).
## Warning: Removed 5 rows containing missing values (geom_point).
Timing of max swe
Monthly max from Nov-Apr
## Warning: Removed 3 rows containing non-finite values (stat_smooth).## Warning: Removed 3 rows containing missing values (geom_point).
Comparison with PRISM data
## Parsed with column specification:
## cols(
## Date = col_date(format = ""),
## ppt = col_double()
## )## Parsed with column specification:
## cols(
## Date = col_date(format = ""),
## tmean = col_double()
## )
freezing days
## Picking joint bandwidth of 20.6
Temperature of hottest average temperature days each water year
Number of freezing days, when average temperature is at or below 0
freezing days correlation to SNOTEL swe max 
##
## Pearson's product-moment correlation
##
## data: snotel_swe_max$swe_mm and snotel_swe_max$ndayfr
## t = 3.1496, df = 9, p-value = 0.01174
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.2196123 0.9230347
## sample estimates:
## cor
## 0.7240951##
## Call:
## lm(formula = snotel_swe_max$swe_mm ~ snotel_swe_max$ndayfr)
##
## Residuals:
## Min 1Q Median 3Q Max
## -110.40 -75.27 -11.91 37.63 187.91
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -231.049 151.128 -1.529 0.1607
## snotel_swe_max$ndayfr 6.491 2.061 3.150 0.0117 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 107.3 on 9 degrees of freedom
## Multiple R-squared: 0.5243, Adjusted R-squared: 0.4715
## F-statistic: 9.92 on 1 and 9 DF, p-value: 0.01174Freezing days relationship to week of water year of max swe
| Data | Median | Mean | Variance | SD |
|---|---|---|---|---|
| SNOTEL SWE (mm) | 45.72 | 90.31 | 13273.55 | 115.21 |
| SNODAS SWE (mm) | 23.00 | 74.12 | 11970.78 | 109.41 |
| PRISM Accumulated Snowfall (mm) | 60.54 | 82.59 | 6924.39 | 83.21 |
##
## Call:
## lm(formula = wc_prism_max_lm$snow_mm ~ wc_prism_max_lm$snotel_swe_mm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -34.480 -12.115 -2.662 8.722 46.558
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 66.91982 15.68070 4.268 0.00209 **
## wc_prism_max_lm$snotel_swe_mm 0.52618 0.05744 9.160 0.00000739 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 26.81 on 9 degrees of freedom
## Multiple R-squared: 0.9031, Adjusted R-squared: 0.8924
## F-statistic: 83.9 on 1 and 9 DF, p-value: 0.000007391##
## Call:
## lm(formula = wc_snodas_max_lm$snow_mm ~ wc_snodas_max_lm$snotel_swe_mm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -66.651 -5.405 11.249 14.802 25.752
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -14.35588 16.60040 -0.865 0.41
## wc_snodas_max_lm$snotel_swe_mm 1.00035 0.06081 16.449 0.0000000505
##
## (Intercept)
## wc_snodas_max_lm$snotel_swe_mm ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 28.38 on 9 degrees of freedom
## Multiple R-squared: 0.9678, Adjusted R-squared: 0.9642
## F-statistic: 270.6 on 1 and 9 DF, p-value: 0.00000005051## [1] 29.37284## Warning: Removed 25 rows containing missing values (geom_point).
##
## Pearson's product-moment correlation
##
## data: wc_prism_snotel$snotel_av_temp and wc_prism_snotel$tmean
## t = 704, df = 4008, p-value < 0.00000000000000022
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9957247 0.9962218
## sample estimates:
## cor
## 0.9959809## Warning: Removed 1 rows containing missing values (geom_path).
## [1] 67.06937##
## Call:
## lm(formula = wc_prism_snotel_wint$prism_snow_mm_accum ~ wc_prism_snotel_wint$snotel_swe_mm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -62.83 -33.97 -10.17 15.80 196.90
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 39.838709 1.522472 26.17
## wc_prism_snotel_wint$snotel_swe_mm 0.581905 0.009901 58.77
## Pr(>|t|)
## (Intercept) <0.0000000000000002 ***
## wc_prism_snotel_wint$snotel_swe_mm <0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 46.55 on 1659 degrees of freedom
## Multiple R-squared: 0.6756, Adjusted R-squared: 0.6754
## F-statistic: 3454 on 1 and 1659 DF, p-value: < 0.00000000000000022## [1] 67.06937
## [1] "Adjusted R2 value is"## [1] 0.8923609precipitation comparison
Beaver Spring
Station ID: 1143
State: Arizona
Network: SNOTEL
County: Apache
Elevation: 9200 ft.
Latitude: 36.33
Longitude: -109.06
HUC: 140802040401
First explore NAs
## # A tibble: 8 x 3
## variable n_miss pct_miss
## <chr> <int> <dbl>
## 1 snow_depth_mm 50 1.28
## 2 air_temp_min_C 30 0.770
## 3 air_temp_av_C 29 0.744
## 4 air_temp_max_C 29 0.744
## 5 Date 0 0
## 6 precip_mm 0 0
## 7 swe_mm 0 0
## 8 station 0 0
Daily basic trends
## Warning: Removed 30 rows containing missing values (geom_path).## Warning: Removed 29 rows containing missing values (geom_path).
Monthly basic trends
Whiskey creek vs Beaver springs
## Warning: Removed 30 rows containing missing values (geom_path).
Seasonal trends for both snotel sites
## Warning: Removed 156 rows containing non-finite values (stat_smooth).## Warning: Removed 156 rows containing missing values (geom_point).## Warning: Removed 223 rows containing non-finite values (stat_smooth).## Warning: Removed 223 rows containing missing values (geom_point).
timing of peak swe
SNOTEL vs SNODAS
Daily data
## Joining, by = "Date"## Warning: Removed 2049 rows containing missing values (geom_path).
## Adding missing grouping variables: `month`
## Warning: Removed 2055 rows containing missing values (geom_path).
## Warning: Removed 8 rows containing non-finite values (stat_smooth).## Warning: Removed 8 rows containing missing values (geom_point).
## Warning: Removed 8 rows containing non-finite values (stat_smooth).
## Warning: Removed 8 rows containing missing values (geom_point).
##
## Call:
## lm(formula = bs_snodas_snotel_wint$snodas_swe_mm ~ bs_snodas_snotel_wint$snotel_swe_mm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -195.405 -11.368 -6.668 8.008 260.892
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 6.668346 1.112497 5.994
## bs_snodas_snotel_wint$snotel_swe_mm 0.960421 0.007601 126.351
## Pr(>|t|)
## (Intercept) 0.00000000242 ***
## bs_snodas_snotel_wint$snotel_swe_mm < 0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 39.07 on 1991 degrees of freedom
## Multiple R-squared: 0.8891, Adjusted R-squared: 0.8891
## F-statistic: 1.596e+04 on 1 and 1991 DF, p-value: < 0.00000000000000022## [1] 39.43443
Comparison with PRISM data
## Parsed with column specification:
## cols(
## Date = col_date(format = ""),
## ppt = col_double()
## )## Parsed with column specification:
## cols(
## Date = col_date(format = ""),
## tmean = col_double()
## )
freezing days
Temperature of hottest average temperature days each water year
Number of freezing days, when average temperature is at or below 0
# raw prism snowfall
bs_prism_snow <- bs_prism_freeze %>%
select(-waterYear) %>%
addWaterYear() %>%
mutate(month = month(Date)) %>%
drop_na(ppt) %>%
group_by(waterYear) %>%
mutate(prism_snow_mm_accum = cumsum(prism_snow_mm)) %>% # accumulate prism precip
mutate(prism_snow_mm_accum = ifelse(month(Date) %in% c(4,5,6,7,8,9,10), 0,
prism_snow_mm_accum))freezing days correlation to SNOTEL swe max 
##
## Pearson's product-moment correlation
##
## data: snotel_swe_max$swe_mm and snotel_swe_max$ndayfr
## t = 2.8648, df = 9, p-value = 0.01863
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.1549281 0.9124708
## sample estimates:
## cor
## 0.6906187##
## Call:
## lm(formula = snotel_swe_max$swe_mm ~ snotel_swe_max$ndayfr)
##
## Residuals:
## Min 1Q Median 3Q Max
## -125.31 -51.99 -11.22 22.42 176.69
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -152.328 139.245 -1.094 0.3024
## snotel_swe_max$ndayfr 5.365 1.873 2.865 0.0186 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 96.15 on 9 degrees of freedom
## Multiple R-squared: 0.477, Adjusted R-squared: 0.4188
## F-statistic: 8.207 on 1 and 9 DF, p-value: 0.01863PRISM Snow accumulation vs SNOTEL SWE
##
## Call:
## lm(formula = bs_prism_snotel_wint$prism_snow_mm_accum ~ bs_prism_snotel_wint$snotel_swe_mm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -109.15 -30.38 -17.05 19.34 248.51
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 29.76708 1.95403 15.23
## bs_prism_snotel_wint$snotel_swe_mm 0.79082 0.01279 61.83
## Pr(>|t|)
## (Intercept) <0.0000000000000002 ***
## bs_prism_snotel_wint$snotel_swe_mm <0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 54.66 on 1659 degrees of freedom
## Multiple R-squared: 0.6974, Adjusted R-squared: 0.6972
## F-statistic: 3823 on 1 and 1659 DF, p-value: < 0.00000000000000022| Data | Median | Mean | Variance | SD |
|---|---|---|---|---|
| SNOTEL SWE (mm) | 76.20 | 101.06 | 11526.41 | 107.36 |
| SNODAS SWE (mm) | 54.00 | 93.46 | 13775.25 | 117.37 |
| PRISM Accumulated Snowfall (mm) | 68.85 | 98.13 | 10146.48 | 100.73 |
## [1] 63.38916##
## Call:
## lm(formula = bs_prism_max_lm$snow_mm ~ bs_prism_max_lm$snotel_swe_mm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -66.321 -37.817 -8.161 33.669 88.401
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 75.8424 36.1926 2.096 0.065593 .
## bs_prism_max_lm$snotel_swe_mm 0.6667 0.1358 4.910 0.000837 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 54.16 on 9 degrees of freedom
## Multiple R-squared: 0.7281, Adjusted R-squared: 0.6979
## F-statistic: 24.1 on 1 and 9 DF, p-value: 0.0008365##
## Call:
## lm(formula = bs_snodas_max_lm$snow_mm ~ bs_snodas_max_lm$snotel_swe_mm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -38.226 -27.618 -6.703 30.295 54.183
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -45.6611 22.5995 -2.02 0.0741 .
## bs_snodas_max_lm$snotel_swe_mm 1.1732 0.0848 13.84 0.000000227 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 33.82 on 9 degrees of freedom
## Multiple R-squared: 0.9551, Adjusted R-squared: 0.9501
## F-statistic: 191.4 on 1 and 9 DF, p-value: 0.0000002272
## [1] 59.23187
## [1] "Adjusted R2 value is"## [1] 0.6979208PRISM snow accumulation using SNOTEL temperature data
precipitation comparison
# combine with snotel data
bs_prism_snotel <- bs_prism_precip %>%
mutate(Date = ymd(Date)) %>%
rename(prism_ppt_mm = ppt) %>%
left_join(bs_clean_metric, by = "Date") %>%
select(Date, prism_ppt_mm, "snotel_ppt_mm" = precip_mm, snow_depth_mm, swe_mm)